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摘 要 : 目前 对 地 震 前 兆 异 常 的 研究 主要 集中 在 “ 热 ” 和 “ 电 ” 等 方面 ， 很 少 涉及 基准 站 的 GPS 数据 。 然 而 ， 已 经 

学 者 证 明 震 中 附近 基准 站 的 GPS 时 间 序 列 坐标 数据 中 也 草 含 着 大 地 震 的 前 兆 信息 。 针 对 2001~2010 年 间 美国 本 土 发 生 
的 具有 代表 性 的 三 个 地 震 进 行 了 研究 ， 首 次 将 Martingale 理论 运用 于 GPS 数据 处 理 ， 提 出 一 种 异常 提取 算法 ， 进 而 对 
地 震 前 后 ， 震 中 附近 多 个 基准 站 的 GPS 数据 进行 分 析 。 实 验 结果 表明 算法 能 够 有 效 的 反映 大 地 震 前 后 GPS 数据 中 异 
常 的 变化 趋势 ， 为 使 用 GPS 数据 对 大 地 震 进 行 预报 提供 了 更 多 可 能 。 
关键 词 : 地 震 ; Martingale 理论 ; 异常 提取 


中 图 分 类 号 : TP391 


doi: 10.3969/j.issn.1001-3695.2017.09.0931 


Detection of anomalies in GPS data before and after earthquakes 
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Abstract: Most research studies on forecasting of earthquakes focus on thermal or electrical aspects. Scholars have proved GPS 


time series also contains precursory information for large earthquakes. For the problems mentioned before, three typical 


earthquakes hit America from 2001 to 2010 are investigated by using GPS data of several GPS reference stations near the 


epicenters. An algorithm based on Martingale theory is proposed to detect the anomalies in GPS data. The experimental results 


show that the proposed algorithm can effectively reflect the change process in GPS data before and after large earthquakes, 


which offers more possibilities for earthquake forecast from GPS data. 
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0 引言 研究 表明 , 基准 站 的 GPS 时 间 序 列 坐 标 数据 经 过 前 期 处 理 ， 
可 以 求 出 每 日 的 坐标 位 移 解 ， 能 反映 中 、 低 频 地 形 形 变 信息 ， 

地 震 是 一 种 常见 的 自然 灾害 , 它 会 对 震 区 人 民 群 众 的 生命 且 具 有 比较 好 的 稳定 性 ,Wang 等 人 中 使 用 Molchan 图 表 法 ， 
和 财产 造成 巨大 的 损害 。 在 过 去 的 数 十 年 中 ， 不 少 学 者 致力 于 ”已 经 证 明 GPS 形变 中 列 含 着 地 震 的 前 兆 信息 。Sagiya 等 人 109 以 
研究 地 震 前 兆 异 常 ， 以 对 地 震 未 来 的 发 展 趋势 进行 预测 ， 涉 及 ”多 个 GPS 基准 站 时 序 坐 标 数据 为 对 象 , 来 研究 地 震 等 其 他 地 壳 
的 指标 主要 包括 “ 热 ” 上 和 “ 电 ” 白 等 方面 。 例 如 ，Kong 等 人 ”变化 。Yeha 等 上 基于 多 个 GPS 基准 站 建立 一 个 GPS 网 络 ， 发 
BHH NOAA 卫星 遥感 射出 长 波 辐射 (outgoing long-wave Ii GPS 时 序 坐 标 数据 在 2013 年 台湾 南投 和 花莲 地 震 前 出 现 了 
radiation, OLR), 探索 汶川 地 震 和 庐山 地 震 前 OLR 数据 的 异常 。 ”一 定 程度 的 扰动 。 文 献 [12] 中 尝试 使 用 一 种 非 线 性 过 滤器 来 提 
Han 等 人 外 统 计 了 2001-2010 年 日 本 Kakioka 站 台 的 超 低 频 段 — 取 GPS 数据 中 可 能 与 地 震 有 关 的 异常 信号 。 赵 等 (3 利用 中 国 大 
E 磁 发 射 数据 , 试图 找 出 数据 中 异常 和 地 震 之 间 的 关联 Nemec ” 陆 几 个 邻近 GPS 基准 站 数据 ,研究 了 日 本 2011 年 发 生 的 9.0 级 
等 人 5 使 用 统计 方法 研究 地 震 前 低频 电磁 波 数 据 ， 发 现在 地 震 th 震 的 震 前 变化 ， 认 为 震 前 运动 速度 变 慢 的 趋势 可 能 是 地 震 的 
前 短 时 间 内 固定 频率 电磁 波 会 出 现 异常 ,Liu 等 人 [研究 了 台湾 前兆。 因此 ,尽管 将 GPS 数据 作为 地 震 前 兆 异常 中 指标 的 研究 
也 区 1993-2004 年 电离 层 总 电子 含量 异常 和 地 震 的 相关 性 。 还 是 相对 较 少 ， 对 其 进行 分 析 ， 寻 找 其 中 蕴含 的 地 震 前 光 异 常 
Kuo 等 人 中 对 日 本 2011 Tohoku-Oki Mw 9.0 级 地 震 前 ， 电 离 层 言 息 具 有 可 行 性 。 然 而 ， 现 有 成 果 主 要 局 限于 地 理学 科 ， 结 果 
总 电子 含量 的 预测 值 和 观察 值 进行 了 对 比 ， 以 研究 地 震 前 电离 。 分 析 主 要 是 基于 人 工 观察 ， 直 接受 到 领域 专家 知识 和 经 验 的 影 
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数据 挖掘 的 目的 就 是 从 海量 的 、 不 完全 的 、 有 噪声 的 数据 
中 提取 出 隐 含 的 未 知 的 但 是 又 潜在 有 用 的 知识 和 信息 的 过 程 。 
当前 ， 已 有 部 分 学 者 将 数据 挖掘 中 异常 检测 的 方法 运用 于 分 析 
地 震 前 不 同 指标 的 变化 。 例 如 , Marzocchi 等 人 04 使 用 贝 叶 斯 估 
值 方法 运用 于 分 析 地 震 数据 ， 文 献 [15] 中 使 用 小 波 变换 用 于 寻 
找 地 震 前 OLR 数据 的 异常 ,文献 [16] 中 使 用 模糊 神经 算法 用 于 
检测 震 前 电磁 异常 ， 文 献 [17] 中 提出 一 种 基于 误差 和 关键 点 的 
地 震 前 兆 观 测 数据 异常 挖掘 算法 。 

本 文 的 贡献 主要 在 于 : 目前 涉及 使 用 基准 站 的 GPS 数据 来 
分 析 震 前 异常 的 研究 还 比较 少 ， 方 法 也 主要 是 基于 人 工 观察 ， 
对 于 异常 的 判定 不 确定 性 较 大 。 因 此 , 本文 使 用 Nevada Geodetic 
Laboratory 提供 的 GPS 数据 , 基于 Martingale 理论 1 与 GPS 数 
据 特 点 ， 提 出 了 一 种 异常 提取 算法 ADM (anomaly detection 
based on martingale theory for GPS data)。 针 对 发 生 在 美国 比较 
具有 代表 性 的 三 个 地 震 ， 运 用 该 算法 ， 对 地 震 前 后 震中 附近 基 
准 站 的 GPS 时 间 序 列 数据 进行 了 分 析 研 究 。 结 果 不 仅 进一步 证 
SKT GPS 数据 中 蕴藏 着 地 震 的 前 兆 信息 , 而 且 为 将 数据 挖掘 算 
法 运用 于 地 理学 科 提 供 了 借鉴 。 


1 ”相关 数据 


本 文 研究 的 对 象 是 从 2001 年 以 来 ， 发 生 在 美国 陆地 上 
度 小 于 60km 的 地 震中 ， 震 级 最 大 的 Ms7.2 级 地 震 ， 以 及 两 
震中 距离 相对 较 近 的 地 震 。 地 震 的 相关 信息 来 自 美国 地 质 调 
体 如 表 


X 


x oW E 


^ E 


局 (USGS, https://earthquake.usgs.gov/earthquakes/), 
1 所 示 。 


del 地 震 相 关 信息 


Date Latitude Longitude Depth Magnitude 
(y-m-d) CN) ( E) (km) (JMA) 
01-02-28 47.149 -122.727 51.80 6.8 
08-02-21 41.144  -114.872 7.90 62 
10-04-04 32.286 -115.295 9.99 72 


使 用 的 GPS 数据 是 GPS 基准 站 站 点 坐标 ， 分 为 东西 向 、 
北 南 向 和 垂直 向 三 个 方向 ， 来 源 于 Nevada Geodetic Laboratory 
的 数据 共享 服务 网 站 (http://geodesy.unr.edu/), 数据 使 用 IGS08 
框架 (1。 综合 考虑 站 点 位 置 以 及 GPS 数据 的 完整 性 , 实验 中 使 
的 GPS 基准 站 的 信息 如 表 2 所 示 。SEAT 站 和 CLOV 站 分 别 
对 应 2001-02-28 与 2008-02-21 这 两 个 发 生 在 不 同时 间 , 震中 位 
相对 较 近 的 地 震 ， 以 分 析 不 同时 期 ， 同 一 地 理 位 置 附近 ， 地 
震 前 后 基准 站 的 GPS 数据 。MEXI、IID2 与 P500 这 三 个 站 点 
对 应 的 是 发 生 在 2010-04-04 的 M7.2 级 地 震 ， 以 研究 对 于 同一 
地 震 ， 震 中 附近 不 同 基准 站 的 GPS 数据 是 否 出 现 相似 的 异常 。 


S 


2 ”研究 方法 
以 Data={di,4d2,…,d1-1} 表示 某 基 准 站 指定 方向 上 的 GPS Dj 
史 数 据 ，di 表示 当前 数据 。 显 然 ， 当 地 质 板块 的 活动 较为 平稳 


时 ， Darta 中 的 样本 分 布 应 该 维持 在 一 个 相对 稳定 的 水 平 ， 样 本 
之 间 具 有 某 种 相似 的 特征 。 因 此 , 结合 文献 [3], 可 以 将 Data od, 
看 成 按照 时 间 顺 序 ， 以 数据 流 形式 到 来 的 ， 来 自 同一 个 聚 类 中 
的 样本 进行 分 析 。 


表 2 GPS 站 点 相关 信息 


Latitude Longitude Height 
Name 
CN) CE) (m) 
SEAT 47.654 -122.309 44.211 
CLOV 40.558 -114.874 2084.911 
MEXI 32.633 -115.476 -22.339 
IID2 32.706 -115.032 9.792 
P500 32.690 -115.300 -20.595 


承 上 所 述 ， 基 于 最 小 化 同一 聚 类 中 各 样本 差异 原则 PO， 设 
C 为 Daa H RKP O, UR (1) 计算 : 


t-l 
1 
25537 d) 


i-l 


那么 di 对 于 C 的 偏 移 程度 可 以 用 式 〈2) 衡量 : 
s, = s(Data,d;) -] d; —C | Q) 
其 中 el 可 以 是 任意 的 相 异 性 度量 函数 , 实验 中 使 用 欧式 距离 。 
考虑 di. d2.....di- 的 分 布 情况 , 可 以 依据 Data rfi d; 近邻 对 象 
的 个 数 来 衡量 其 异常 程度 。 因 此 ，ADM 算法 中 用 式 (3) 来 衡 
量 d, 所 对 应 的 当前 指定 方向 上 GPS 数据 的 平稳 程度 : 


#{j|sj; >s +l sj =s} 
= 
t 


lim 


G) 


在 式 G) F, 6 是 (0 内 的 随机 数 ，#{} 是 一 个 函数 ， 返 
可 满足 给 定 条 件 的 样本 的 数量 ， JR QD 计 
5. 

从 式 (3) 可 以 看 出 ，pi e(o1] H pi RK, 表明 d, 越 符合 之 
前 样本 的 分 布 情况 ， 相 应 的 当天 数据 出 现 异常 的 可 能 性 越 小 。 
然而 ，GPS 数据 中 不 可 避免 地 会 存在 噪声 (闪烁 噪声 以 及 随机 
漫步 噪声 等 )， 某 天 相对 较 小 的 产值 并 不 足以 说 明 当 前 GPS 数 
据 整 体 出 现 了 异常 ， 需 要 通过 使 用 Martingale 理论 0， 综合 
虑 整个 时 间 序 列 数据 的 变化 趋势 。 

Martingale 理论 0 起 源 于 赌博 业 和 概率 论 ， 是 最 早 的 金融 


j=1,2,..…,t，54 采 


资产 定价 模型 之 一 ， 现 在 已 经 被 广泛 运用 于 决策 优化 ， 投 资 分 
析 等 领域 。 为 了 真实 反映 GPS 数据 持续 变化 的 情况 ， 结 合 


Martingale 理论 与 GPS 数据 特点 ， 所 需要 分 析 的 每 个 数据 di 所 
对 应 的 随机 化 的 Martingale 值 ( 以 下 简称 MM 值 ) 以 式 (4) 计 
算 : 


t 
n 7 [ [Pma f^ (4) 
i-l 


根据 文献 [21]， 取 ==0.82 并 初始 化 Mo =1。 从 式 (4) P M 
值 的 计算 过 程 可 以 看 出 , 越 大 的 M 值 表示 之 前 一 定时 间 内 GPS 
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数据 的 异常 越 明显 。 


聚 类 算法 中 ， 初 始 聚 类 中 心 的 选择 会 对 聚 类 结果 造成 较 大 
的 影响 。 类 似 的 ，Dara 中 首 个 样本 di 位 置 的 选取 也 会 对 异常 的 
提取 造成 较 大 的 影响 。 因 此 ， 算 法 中 设置 了 一 个 初始 聚 类 样本 


num 


在 初 始 化 C2—— V ds 以 及 


k=1 


个 数 参 数 nm, 


si - Mi - Ci 12. num E, M, dama 开始 计算 M 值 。 直 观 地 说 ， 
就 是 将 所 分 析 数据 前 mum 天 的 平均 数据 作为 初始 中 心 ， 并 从 
num 天 后 开始 分 析 以 减少 结果 的 不 确定 性 。 
此 外 ,由 于 大 地 震 发 生前 后 地 过 运动 相对 剧烈 , 表现 为 GPS 
数据 有 可 能 在 短 时 间 内 出 现 较 大 波动 。 结 合式 CAD, MEI M f 
会 迅速 不 断 增 大 到 不 可 控 的 程度 。 为 了 避免 这 种 情况 的 发 生 ， 
算法 中 设置 了 另 一 个 停止 参数 六 。 也 就 是 说 , 如 果 Mi >h, 那么 
从 di 开始 重新 开始 算法 。 参 数 mum 与 有 对 算法 效果 的 影响 将 


在 实验 部 分 进行 进一步 分 析 与 验证 。 
综 上 所 述 ，ADM 算法 的 具体 步骤 如 下 : 
算法 名 称 : ADM 
Wü A 某 基准 站 指定 方向 上 的 GPS 数据 Data 
停止 参数 h ， 
初始 聚 类 样本 数 num 
输 ”出 : Data 中 每 个 样本 对 应 的 M 值 


num 
Step.1 初始 化 样本 中 心 : Crum =— X di 
i-l 
偏 移 程 度 ， sp-|d;-C|i 12... num 
M fli: Mi =1i=1,2,.…,num 
当前 处 理 样本 序号 : konumel 
循环 停止 条 件 : 
当前 批 次 样本 开始 序号 :front=1 
Il 如 果 检测 到 Mj >h, IWAMA dia 开始 
/ 一 个 新 的 批 次 〈《 即 算法 重新 启动 ) 
Step.2 while k<|Datd && ! isStop 
2. fF] dk, Cki 与 式 (2)， 计 算 sk 
2.2 [EJ] S front» 5 front-1«-5k 与 式 (32, 


isStop- false 


计算 Pk 


2.3 使 用 Mk-1 与 式 (4)， 计 算 MK 
2.4 计算 d fronts d front+ HE dy 的 聚 类 中 心 


大 
Ck = > d 
k k — front k 


i= front 
2.5 使 Ck ~ d front» d frontal je dk HX (2)， 更 新 
S front» 8 fronts: 3K 
2.6 k++ 


2.7 if My4»h 


isStop = true 


end if 


end while 


Step.3 if isStop-true && k - num &|Datd 


3.] isStop- false 


3.2 Ckinum = 


i k+num 
di 
num 
i=k 


33 Mi =1l,i=k,k+1,...,k +num 


3.4 front=k 


3.5 k=k+num+1 
3.6 转 到 Step.2 /W 开 始 一 个 新 的 批 次 
end if 


Step4. $H 


3 ”结果 与 分 析 


H Mj,i -12,...|Datd 


a) SEAT 站 和 CLOV 站 的 原始 GPS 时 间 序 列 坐 标 数 据 以 
及 对 应 的 M 值 如 图 1 和 2 所 示 。 两 个 基准 站 分 别 对 应 两 个 发 生 


在 不 同时 间 、 震 中 位 


num 设 置 为 45,h 设置 为 


10 月 开始 ,一直 延续 到 
值 对 结果 造成 的 影响 ， 医 
的 观察 GPS i3) 


相对 较 近 的 地 震 。ADM 算法 中 ， 参 数 
1000。 研 究 的 数据 从 发 生地 震 上 年 的 
也 震 后 的 一 个 月 。 为 了 降低 2 的 随机 取 
FP 使 用 10 次 实验 的 平均 值 。 出 于 更 好 
居 变 化 趋势 的 目的 ， 原 始 值 中 只 使 用 相应 GPS 


数值 的 小 数 部 分 而 忽略 相同 的 整数 部 分 ， 黄 色 竖 线 表 示 地 震 发 


生 时 间 。 


从 图 1 (d) 与 
值 几 乎 没有 任何 变化 , 说 明 
两 个 基准 站 东西 向 和 北 南 
地 震 后 一 周 左右 迅速 达到 


图 2 (d) 可 以 看 出 ， 前 期 各 个 方向 上 的 M 


c 


GPS 数据 比较 稳定 ,临近 地 震 前 后 ， 


向 的 W 值 均 出 现 了 显著 变化 , 并 且 在 
E, 这 进一步 证 实 了 GPS 数据 中 药 


含 着 和 地 震 有 关 的 信息 。 f 
动 ， 在 地 震 后 才 达到 峰值 


为 了 反映 GPS 


直 得 注意 的 是 ，M 值 在 地 震 前 出 现 波 
的 原因 可 能 在 于 : 于 地 震 的 生成 机 


jj, GPS 数据 的 异常 可 能 
2803] 


出 现在 地 震 前 很 短 的 一 段 时 间 内 ， 而 
化 趋势 以 及 减少 噪声 对 结果 造成 的 影 


出 现 显著 的 变化 ， 
的 表现 为 地 震 后 短 时 间 


响 ， 算 法 需要 在 累积 一 定数 量 的 异常 样本 之 后 ，M 值 才 会 相应 


天 


内 M 值 达 到 峰值 。 虽 然 可 以 通过 将 参数 
设置 为 较 小 的 数值 以 将 M 值 峰值 出 现 的 时 间 提 前 ,但 是 这 有 


可 能 导致 出 现 误 报 , 本 文 将 在 实验 co 中 验证 参数 对 算法 效果 


的 影响 。 


由 于 两 个 地 震 震 中 位 置 相对 较 近 ， 成 因 相似 ， 因 此 所 提取 


的 东西 向 和 北 南 向 的 W 值 变化 规律 也 比较 类 似 , 这 从 侧面 上 说 


明 本 文 提出 的 ADM 算法 


的 可 靠 性 。 此 外 ， 这 两 个 方向 上 的 原 


台 GPS 数据 图 
法 提取 GPS 数 


据 异 常 的 


误差 要 远 远 大 了 


FR 


像 并 无 出 现 明显 的 异常 趋势 ,这 证 实 了 ADM 算 
了 效 性 。 因 为 垂直 向 GPS 数值 的 测量 
向 和 北 南 向 名 ， 所 以 对 比 东西 向 和 北 南 向 
的 W 值 , 两 个 基准 站 垂直 向 M 值 波动 都 不 明显 , 并 没 
地 震 明 显 相 关 的 规律 。 
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图 2 


b) 2010-04-04 地 震 震 中 邻近 MEXI 站 、IID2 站 和 P500 站 ” 别 单独 列 出 对 应 三 个 方向 上 的 M 值 。 
的 GPS 数据 原始 值 以 及 对 应 的 M 值 分 别 见 图 3、 图 4 与 图 5， 对 比 图 3(c) (qd) 可 以 看 出 , 由 于 地 震 后 基准 站 特定 方向 上 的 
原始 值 中 同样 只 显示 相应 的 小 数 部 分 ， 黄 色 竖 线 表示 地 震 发 生 ”GPS 数据 出 现 的 “断崖 ” 式 的 变化 ， 相 应 的 M 值 也 在 之 后 对 应 
时 间 。 实 验 的 目的 是 研究 对 于 同一 地 震 ， 震 中 附近 不 同 基 准 站 ”的 一 段 时 间 内 达到 峰值 。 类 似 的 情况 也 出 现在 图 4(a)(b)， 图 
GPS 数据 的 异常 。ADM 算法 中 ， 参 数 mm 设置 为 45,h 设置 为 。 5(c)(d) 等 。 波峰 的 延 后 出 现 , 一 方面 说 明 本 文 提出 的 ADM 算法 
1000。 由 于 此 次 地 震 后 短 时 间 内 ,各 个 基准 站 特定 方向 上 的 GPS ”能 够 真实 反映 GPS 数据 中 的 异常 , 男 一 方面 也 说 明了 实验 1 中 
数据 会 突然 出 现 “ 断 崖 ” 式 的 变化 ， 因 此 研究 的 数据 从 发 生地 ”导致 M 值 出 现 峰值 的 GPS 数据 异常 可 能 是 在 震 前 出 现 的 。 
震 上 年 的 10 月 开始 , 一 直 延 续 到 地 震 后 的 两 个 月 。 为 了 更 清晰 1(9) 与 图 1(b) 中 , MEXI 站 东西 向 没有 出 现 类 似 情况 的 原因 在 于 : 
的 展示 地 震 前 后 , 三 个 方向 GPS 数据 的 异常 趋势 , 图 3~5 中 分 GPS 数据 在 地 震 前 就 出 现 了 潜在 的 剧烈 变化 , 导致 M 值 超过 了 
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预先 设置 的 停止 参数 =1000, 于 是 算法 重新 启动 


分 


大 ， 


的 M 
GPS 数据 就 有 较 大 的 异常 发 生 , 从 侧 


拟 合 了 这 部 


“断崖 ”数据 。 


可 能 是 由 于 地 震 成 因 不 同 并 且 2010-04-04 地 震 的 震级 较 
丸 此 三 个 基准 站 东西 向 的 M 值 , MEXI 站 与 ID2 站 北 南 向 
也 震 前 一 段 时 间 就 达到 过 峰值 , 说 明 在 地 震 发 生前 ” ”可 能 与 基准 站 的 具体 位 置 有 关 。 与 实验 
下 证 明了 GPS 形变 中 蕴含 E 


值 都 在 


着 地 震 的 前 兆 信息 , 尽管 此 时 原始 GPS 数据 
常 趋势 。 这 为 使 用 本 文 提 出 的 ADM 算法 ,利用 GPS 数据 对 大 


又 出 现 平缓 的 原因 


站 重 


EF: 算法 重新 启动 ， 需 要 拟 合 mum 天 的 数 
据 后 再 判定 异常 ， 这 段 时 间 内 M. 值 始终 为 初始 值 1。 与 MEXI 


站 与 IID2 站 类 似 ，P500 站 北 南 向 的 NM 值 在 地 震 前 也 出 现 了 一 


定 程度 的 变化 ， 只 不 过 变化 程度 要 远 低 于 上 述 两 个 基准 站 ， 这 


旭 像 并 无 明显 的 异 ”垂直 向 的 NM 值 比较 明显 


地 震 进 行 预报 提供 了 可 能 。 各 个 基准 站 第 
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(b) MEXI 站 东西 向 M 值 (2010-04-04 地 震 ) 
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(e) MEXI 站 垂直 向 原始 值 (2010-04-04 地 震 ) 
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(e)IID2 站 垂直 向 原始 值 (2010-04-04 地 震 ) 


图 4 


AH M 值 的 变化 并 没有 呈现 一 致 的 规律 ， 
也 在 地 震 前 后 出 现 了 ; 


4 有 IID2 站 
ZJ, MEXI xl; 


NI 


直 向 的 波动 可 能 来 源 于 数据 “断崖 ” 式 的 变化 ,而 P500 站 垂直 
向 M 值 的 峰值 出 现在 地 震 后 一 段 时 间 内 。 
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图 5 
从 基准 站 分 布 的 位 置 上 看 ， 图 3~5 中 可 以 看 出 ， 由 于 三 个 。” 小 的 2001-02-28 Hij (Ms6.8) 与 2008-02-21 地 震 (Ms6.2) 各 


基准 站 都 位 于 震中 附近 ， 因 此 在 测量 误差 相对 较 小 的 东西 向 和  ” 自 对 应 基准 站 的 M 值 则 最 早 在 地 震 前 一 个 月 左右 开始 变化 , 并 
北 南 向 上 , 各 个 基准 站 M 值 第 一 次 出 现 波峰 的 时 间 并 不 存在 显 。 在 地 震 后 较 短 的 一 段 时 间 内 达到 波峰 。 这 说 明 地 震 震 级 越 大 ， 
著 差 异 。 但 是 ， 三 个 基准 站 中 距离 震中 最 近 的 MEXI 站 两 个 方 之 前 需要 积 揭 的 能 量 越 多 ,GPS 数据 潜在 的 异常 可 能 更 早 发 生 ， 
向 上 波峰 第 一 次 出 现 的 时 间 会 早 于 IID2 站 与 P500 站 。 这 说 明 KEM 值 相应 越 早出 现 波动 。 此 外 ，2001-02-28 地 震 (Ms6.8) 


距离 震中 越 近 的 基准 站 ，MM 值 可 能 越 早 出 现 波 峰 ， 即 GPS 数 M 值 波 动 出 现 的 时 间 要 述 于 2008-02-21 地 震 (Ms6.2) 的 原因 
据 潜 在 异常 越 早 发 生 。 可 能 在 于 前 者 的 震源 深度 (51.80km) 要 远大 于 后 者 (7.90km )， 
结合 图 1~5， 同 样 从 测量 误差 相对 较 小 的 东西 向 和 北 南 向 因此 地 震 发 生前 期 GPS 数据 潜在 变化 相对 不 明显 。 


E, 进一步 分 析 地 震 震 级 、 发 生 时 间 与 M 值 之 间 的 关系 。 本 文 c) 为 了 评价 参数 有 对 算法 效果 的 影响 ， 本 文 在 MEXI 站 北 
研究 的 三 个 地 震中 ， 震 级 最 大 的 2010-04-04 地 震 (Ms7.2) 所 WEE, 分 别 将 4 设置 为 250,500,1000,2000 进行 了 测试 ， 参数 
对 应 的 三 个 基准 站 的 M 值 第 一 次 波峰 出 现 的 时 间 均 在 地 震 发 。 num 设置 为 45。 实 验 结果 如 图 6 所 示 。 

生前 三 个 月 左右 ， 变 化 一 直 持 续 到 地 震 发 生 时 。 而 震级 相对 较 
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图 6 hh 值 对 算法 效果 的 影响 (MEXI 站 北 南 向 ) 


从 图 6 中 可 以 看 出 ， 当 GPS 数据 中 出 现 异常 时 ，M 值 将 。 会 有 较 大 的 差异 ， 这 有 可 能 导致 波峰 出 现 沾 后 。 由 于 过 小 的 4 
以 指数 级 别 迅速 增长 。 显然 , 当 num 相 同时 , 通过 将 参数 4 设置 值 有 可 能 会 导致 过 多 的 误 报 , 因此 h=1000 是 一 个 相对 合理 的 取 
为 较 小 的 数值 ， 可 以 将 首次 M 值 峰值 出 现 的 时 间 提 前 。 对 比 — 值 。 
500,1000,2000， 将 设置 为 250 时 ， 第 二 次 波峰 出 现 的 时 间 并 dO 为 了 评价 参数 num 对 算法 效果 的 影响 ， 本 文 同样 在 
没有 提前 的 原因 在 于 : 当 第 一 次 检测 到 M r, 算法 需要 重新  MEXI 站 北 南 向 上 , 分 别 将 num 设置 为 25,35,45,55 进行 了 测试 ， 
启动 。 此 时 ,不 同 4 值 重新 开始 的 位 置 不 同 ， 在 GPS 数据 中 已 ”参数 疡 设置 为 1000。 实 验 结果 如 图 7 所 示 。 
经 出 现 异 常 的 前 提 下 ， 再 次 获得 的 新 批 次 样本 的 初始 聚 类 中 心 


201805.00377v1 


chinaXiv 


录用 稿 


2000 


ChinaXiv& fER 1 
李 h, 等 : 地 震 前 后 基准 站 VERAT 法 


1500 


1000 


M Value 


500 


2009.75 2009.87 2009.99 


R] 


显然 , 24 GPS 数据 稳定 时 , 只 要 num 足够 大 , 受 其 影响 的 
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Date 


7 _num 值 对 算法 效果 的 影响 IMEXI 站 北 南 向 ) 


2010.23. 


2010.37 


2010.49 


着 


地 震 的 前 兆 信息 ， 


也 为 使 用 GPS 数据 对 地 震 进 行 预报 提供 了 


不 同 wm 值 所 对 应 的 初始 聚 类 中 心 C 差别 应 该 不 大 。 因 此 ， 从 可 能 。 然 而 ， 基 于 目前 的 数据 量 和 研究 水 平 ， 探 究 使 用 基准 站 
图 7 中 可 以 看 出 ，num 值 的 不 同 对 M 值 第 一 次 波峰 出 现 的 时 间 GPS 数据 对 大 地 震 进 行 预报 将 是 一 项 长 期 工作 。 本 文 提出 的 
几乎 没有 影响 。 然 而 ， 与 实验 (GO 类 似 ， 当 算法 重新 启动 时 ， ADM 算法 所 提取 的 震中 附近 基准 站 GPS 数据 所 对 应 的 W 值 的 
随 着 GPS 数据 中 潜在 异常 的 出 现 , 不 同 num 值 所 对 应 的 初始 聚 ” 异常 已 经 揭示 了 一 些 初 步 规律 ， 但 还 需要 在 更 多 的 地 震 实例 上 
类 中 心 C 有 所 差别 , KE M 值 第 二 次 出 现 波峰 的 时 间 存 在 一 定 。 ”进一步 验证 M 值 的 变化 趋势 与 地 震 各 参数 震级、 位 置 、 发 生 
的 差异 。 虽然 实验 (4) 中， 将 mm 设置 为 55 时 M 值 出 现 波峰 ”时 间 等 ) 之 间 的 关系 并 对 其 进行 统计 可 靠 性 分 析 ， 这 也 是 本 文 
的 时 间 最 早 ， 但 是 考虑 到 将 来 实际 应 用 中 ， 过 大 的 num 值 有 可 下 一 步 的 研究 方向 。 

Lb d; 2H yr yz | 如 A VS He ` N ut 

uin M 值 两 次 波峰 出 现 的 间隔 过 长 ， 因 此 本 文中 将 num 设 置 参考 文献 ; 

e) 为 了 比较 本 文 提 出 的 ADM 算法 与 目前 常用 的 利用 领域 [1] Saradjiaa M R, Alhoondzadeh M. Thermal anomalies detection before 
ERINRAR Se AL HR BE JE 055] 7I 3C), ARX E SEAT strong earthquakes (M»6, 0) using interquartile, wavelet and Kalman filter 
站 北 南 向 的 数据 上 进行 了 实验 , 数据 异常 判定 采用 ko 准则 。 在 method [J]. Natural Hazards and Earth System Sciences, 2011, 11 (1): 1099- 
ko 准则 中 ， 以 x 表示 观测 结果 的 平均 值 ，o 表示 观测 结果 的 标 1108. 

准 差 , 如 果 数 据 超出 了 [x-ko,x+ko] 的 范围 , 则 认为 有 异常 发 生 [2] Honkura Y, Oshiman N, Matsushima M, et al. Rapid changes in the electrical 

实验 中 依据 文献 [22] 将 大 值 设 定 为 2， 结 果 如 图 SPI. state of the 1999 izmit earthquake rupture zone [J]. Nature Communications, 
2013, 4 (2116): 1-10. 
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